home *** CD-ROM | disk | FTP | other *** search
/ MacWorld 1999 November / Macworld (1999-11).dmg / Updaters / WhiteCap 3.0.4 / WhiteCap Source.sit / WhiteCap Source / Common / General Tools / XFloatList.cpp < prev    next >
C/C++ Source or Header  |  1999-07-13  |  6KB  |  262 lines
  1. #include "XFloatList.h"
  2.  
  3. #include "Hashtable.h"
  4. #include "XLongList.h"
  5. #include <stdlib.h>
  6. #include <math.h>
  7.  
  8. #define _MIN( a, b ) (( (a) > (b) ) ? (b) : (a) )
  9. #define _MAX( a, b ) (( (a) > (b) ) ? (a) : (b) )
  10.  
  11.  
  12.  
  13. XFloatList::XFloatList( ListOrderingT inOrdering ) :
  14.     mList( inOrdering ) {
  15.  
  16.     if ( inOrdering == cSortLowToHigh || inOrdering == cSortHighToLow )
  17.         mList.SetCompFcn( sFloatComparitor, inOrdering == cSortLowToHigh );
  18. }
  19.  
  20.     
  21.     
  22.  
  23.  
  24.  
  25. void XFloatList::FindMeans( long inNumMeans, float outMeans[], float inSigmaScale ) const {
  26.     long start, end, m, i, n = mList.Count();
  27.     float* srce = (float*) mList.getCStr(), *smoothed = new float[ n ], *temp = NULL;
  28.     float sigma = 0.1 + inSigmaScale * ( (float) ( n  / inNumMeans ) ); 
  29.     float left, cen, right, sum;
  30.     
  31.     // We need all out numbers in order, we must sort them if they're no sorted
  32.     if ( mList.mOrdering != cSortHighToLow ) {
  33.         
  34.         // Copy all the array elements for sorting...
  35.         temp = new float[ n ];
  36.         for ( i = 0; i < n; i++ )
  37.             temp[ i ] = srce[ i ];
  38.         
  39.         // Sort all the floats from this to temp
  40.         qsort( temp, n, 4, sQSFloatComparitor );
  41.         srce = temp;
  42.     }
  43.     
  44.     // Smooth all the values (they're already sorted)
  45.     GaussSmooth( sigma, n, srce, smoothed );
  46.     
  47.     // Compute the discrete 1st derivative
  48.     for ( i = 0; i < n - 1; i++ )
  49.         smoothed[ i ] = fabs( smoothed[ i ] - smoothed[ i + 1 ] );
  50.     
  51.     // We want to find the top <inNumMeans> local max of the 1st deravative
  52.     Hashtable sepCandidates;
  53.     cen = smoothed[ 0 ];
  54.     right = smoothed[ 1 ];
  55.     for ( i = 1; i < n - 2; i++ ) {
  56.     
  57.         left = cen;
  58.         cen = right;
  59.         right = smoothed[ i+1 ];
  60.         
  61.         // If this a local max.  (Note: this could/should be improved for neighbors that are equal)
  62.         if ( ( cen > left && cen >= right ) ) {
  63.             sepCandidates.Put( i, *((void**) &cen) );
  64.         }
  65.     }
  66.     
  67.     // Pick out the 1st derative peaks, then dealloc what we no longer need
  68.     XPtrList    rank;
  69.     sepCandidates.Rank( rank, sQSFloatComparitor, inNumMeans - 1 );
  70.     delete []smoothed;
  71.  
  72.     XLongList    quintiles( cSortLowToHigh );
  73.     for ( i = 1; i < inNumMeans; i++ )
  74.         quintiles.Add( (long) rank.Fetch( i ) );
  75.     quintiles.Add( n );
  76.     
  77.     // The means are the averages of the initial (sorted) data, divided by the serparating ranks
  78.     start = 0;
  79.     for ( m = 1; m <= inNumMeans; m++ ) {
  80.         end = quintiles.Fetch( m );
  81.         for ( sum = 0, i = start; i < end; i++ )
  82.             sum += srce[ i ];
  83.         outMeans[ m - 1 ] = sum / ( (float) ( end - start ) );
  84.         start = end + 1;
  85.     }
  86.     
  87.     // Cleanup
  88.     if ( temp )
  89.         delete []temp;
  90. }
  91.     
  92.     
  93.  
  94.  
  95.     
  96. #define _safeSmooth( xx )                                    \
  97.         sum = 0;                                            \
  98.         factor = 1;                                            \
  99.         for ( i = - maskCtr; i <= maskCtr; i++ ) {            \
  100.             rx = xx + i;                                    \
  101.             if ( rx < 0 || rx >= inN )                         \
  102.                 factor -= mask[ i + maskCtr ];                \
  103.             else                                            \
  104.                 sum += mask[ i + maskCtr ] * inSrce[ rx ];    \
  105.         }                                                    \
  106.         inDest[ xx ] = sum / factor;
  107.                 
  108.  
  109.  
  110.  
  111. void XFloatList::GaussSmooth( float inSigma ) {
  112.  
  113.     GaussSmooth( inSigma, mList.Count(), (float*) mList.getCStr() );
  114.     
  115. }
  116.  
  117.  
  118. void XFloatList::GaussSmooth( float inSigma, long inN, float inSrceDest[] ) {
  119.     float* temp = new float[ inN ];
  120.     
  121.     GaussSmooth( inSigma, inN, inSrceDest, temp );
  122.     
  123.     for ( long i = 0; i < inN; i++ )
  124.         inSrceDest[ i ] = temp[ i ];
  125. }
  126.  
  127.  
  128. void XFloatList::GaussSmooth( float inSigma, long inN, float inSrce[], float inDest[] ) {
  129.     int        maskSize        = _MAX( 8.0 * inSigma, 4.0 );
  130.  
  131.     // Make sure the mask size is odd (so we have a 'center')
  132.     if ( maskSize % 2 == 0 ) 
  133.         maskSize++;
  134.     int        i, x, rx, xEnd;
  135.     int        maskCtr            = maskSize / 2;
  136.     float*    mask            = new float[ maskSize ];
  137.     float    sqrt2PiSigma    = sqrt( 2.0 * 3.14159 ) * inSigma;
  138.     float    expon, sum, factor, *base;
  139.  
  140.     // Generate a normalizes gaussian mask out to 4*sigma on each side
  141.     // We have to adjust the center weight so that when sigma falls below 1.0ish so that the sample it 
  142.     // takes at x = 0 doesn't represent the val of a gaussian for -.5 to .5 (when sigma is small
  143.     // a gaussian gets very 'high' at x=0 vs. x=+/- .5 when sigma is small
  144.     sum = 0;
  145.     for ( x = - maskCtr; x <= maskCtr; x++ ) {
  146.         expon = -0.5 * ((float) (x * x)) / ( inSigma * inSigma );
  147.         mask[ x + maskCtr ] = exp( expon ) / sqrt2PiSigma;
  148.         if ( x != 0 )
  149.             sum += mask[ x + maskCtr ];
  150.     }
  151.     
  152.     // This forces normalized weights.
  153.     mask[ maskCtr ] = 1.0 - sum;
  154.         
  155.     // Smooth the lefthand side of the sequence
  156.     xEnd = _MIN( maskCtr, inN );
  157.     for ( x = 0; x < xEnd; x++ ) {
  158.         _safeSmooth( x )
  159.     }
  160.     
  161.     // Smooth the center portion the sequence
  162.     xEnd = inN - maskCtr;
  163.     base = inSrce;
  164.     for ( x = maskCtr; x < xEnd; x++ ) {
  165.         
  166.         // For each pixel, apply the mask to its neighbors for a final value of that pixel
  167.         sum = 0;
  168.         for ( i = 0; i < maskSize; i++ ) 
  169.             sum += mask[ i ] * base[ i ];
  170.         
  171.         base++;
  172.         inDest[ x ] = sum;
  173.     }
  174.     
  175.     // Smooth the righthand side of the sequence
  176.     for ( x = _MAX( maskCtr, inN - maskCtr ); x < inN; x++ ) {
  177.         _safeSmooth( x )
  178.     }
  179.  
  180.     
  181.     // Cleanup
  182.     delete []mask;
  183. }
  184.     
  185.         
  186.  
  187. void XFloatList::Rank( XLongList& outRank, long inNumToRank ) const {
  188.     long *p, *temp;
  189.     float* srce;
  190.     long i, n = Count();
  191.  
  192.     outRank.RemoveAll();
  193.     
  194.     if ( inNumToRank < 0 )
  195.         inNumToRank = n;
  196.     inNumToRank = _MIN( inNumToRank, n );
  197.     
  198.     // Handle trivial cases of this lis already being sorted
  199.     if ( mList.mOrdering == cSortLowToHigh ) {
  200.         for ( i = 0; i < inNumToRank; i-- ) 
  201.             outRank.Add( n - i ); }
  202.             
  203.     // Duh... still sorted...
  204.     else if ( mList.mOrdering == cSortHighToLow ) {
  205.         for ( i = 1; i <= inNumToRank; i++ ) 
  206.             outRank.Add( i ); }    
  207.             
  208.     else {
  209.         temp = new long[ 2 * n ];
  210.         srce = (float*) mList.getCStr();
  211.         
  212.         // To rank, we must sort, with a tag on each element saying where it came from
  213.         p = temp;
  214.         for ( i = 1; i <= n; i++ ) {
  215.             *((float*) p) = *srce;
  216.             p++;  srce++;
  217.             *p = i;
  218.             p++;
  219.         }
  220.         
  221.         // Sort the floats
  222.         qsort( temp, n, 8, sQSFloatComparitor );
  223.         
  224.         // Put the sorted results in the destination
  225.         p = temp + 1;
  226.         for ( i = 0; i < inNumToRank; i++ ) {
  227.             outRank.Add( *p );
  228.             p += 2;
  229.         }
  230.         
  231.         // Cleanup
  232.         delete []temp;
  233.     }
  234. }
  235.  
  236.  
  237.  
  238. int XFloatList::sQSFloatComparitor( const void* inA, const void* inB ) {
  239.     float diff =  *((float*) inB) - *((float*) inA);
  240.     if ( diff > 0.0 )
  241.         return 1;
  242.     else if ( diff < 0.0 )
  243.         return -1;
  244.     else
  245.         return 0;
  246. }
  247.  
  248.         
  249.  
  250. int XFloatList::sFloatComparitor( const void* inA, const void* inB ) {
  251.     float diff =  *((float*) &inB) - *((float*) &inA);
  252.     if ( diff > 0.0 )
  253.         return 1;
  254.     else if ( diff < 0.0 )
  255.         return -1;
  256.     else
  257.         return 0;
  258. }
  259.  
  260.  
  261.  
  262.